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ABSTRACT 

We present equilibrium models of relativistic magnetised, infinite, axisymmetric jets 
with rotation propagating through an homogeneous, unmagnetised ambient medium 
at rest. The jet models are characterised by six functions defining the radial profiles 
of density, pressure, and the toroidal and axial components of velocity and magnetic 
field. Fixing the ambient pressure and the jet rest-mass density and axial components 
of the flow velocity and magnetic field, we analyze the influence of the toroidal mag¬ 
netic field and several rotation laws on the structure of the equilibrium models. Our 
approach excludes by construction the analysis of the self-consistently magnetically 
launched jet models or the force-free equilibrium solutions. Several forbidden regions 
in the magnetic pitch angle/magnetization plane are found where models of the class 
considered in our study could not be settled. These forbidden regions are associated 
with the existence of maximum axial and toroidal magnetic field components compat¬ 
ible with the prescribed equilibrium condition at the jet surface, and/or an excess of 
centrifugal force producing gaps with negative pressures in the jet. The present study 
can be easily extended to jet models with different transversal profiles and magnetic 
field configurations. 

In the last part of the paper, we test the ability of our RMHD code to maintain 
steady equilibrium models of axisymmetric RMHD jets in one and two spatial dimen¬ 
sions. The one dimensional numerical simulations serve also as a consistency proof of 
the fidelity of the analytical steady solutions discussed in the first part of the paper. 
The present study allows us to build initial equilibrium jet models with selected prop¬ 
erties for dynamical (and emission) simulations of magnetised relativistic jets with 
rotation. 

Key words: Galaxies: jets; Physical Data and Processes: magnetic fields, magneto¬ 
hydrodynamics; Methods: analytical, numerical 


1 INTRODUCTION 

Understanding the physics behind the relativistic jets 
emanating from radio-loud AGN (Boettcher, Harris & 
Krawczynski 2012) is a prominent science case that is 
pushing continuously the limits of observational capabili¬ 
ties through the whole electromagnetic spectrum, from ra¬ 
dio to gamma-rays. The Event Horizon Telescope (Clery 
2012), able to achieve 20 /xas resolution at submilimeter 
wavelengths, aimed to capture General Relativistic signa¬ 
tures near the horizon of the central black hole where these 
jets are generated from, is a good example of the continuous 
improvement in observational skills. As important as these 
instrumental/technical advances are high-resolution numer¬ 
ical simulations that, as a sort of virtual laboratory, try to 
connect both theoretical models and more and more precise 
observations. 


For parsec scale jets, the strategy of combining analyt¬ 
ical steady jet solutions (acting as initial conditions) with 
dynamical and emission simulations has proven already its 
success in the relativistic, non-magnetised regime (Gomez et 
al. 1997; Komissarov & Falle 1997; Agudo et al. 2001; Aloy 
et al. 2003; Mimica et al. 2009; Fromm et al. 2012, Fromm 
et al. 2015, in prep.) giving plausible interpretations to most 
of the phenomenology found in these objects. 

Relativistic magnetohydrodynamic simulations have 
concentrated on the morphological characterisation of large 
scale (i.e., kiloparsec scale) magnetised jets (Komissarov 
1999b; Leismann et al. 2005; Keppens et al. 2008; Mignone 
et al. 2010) and the formation of jets from magnetohydro- 
dynamical mechanisms (see, e.g., McKinney & Blandford 
2009; Tchekhovskoy, Narayan & McKinney 2011; McKin¬ 
ney, Tchekhovskoy & Blandford 2012; Porth 2013). Porth et 
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al. ( 2011 ) computed the synchrotron radiation (in the mil¬ 
limeter and submillimeter range) from RMHD jets through 
the acceleration region allowing them to address important 
issues as the frequency-dependent core-shift effect, or the 
signatures of large scale helical helds in the synchrotron radi¬ 
ation. Previously, Roca-Sogorb et ah (2008, 2009) presented 
exploratory results on the effects of helical magnetic fields 
on the dynamics of non-rotating relativistic compact jets 
and their potential observable imprints, also using RMHD 
simulations. Nakamura, Garofalo & Meier (2010) and Naka¬ 
mura k, Meier (2014) interpreted the proper motions of 
components ejected from the HST-1 complex from radio to 
optical wavelengths as magnetohydrodynamic shock fronts 
making use of one-dimensional RMHD simulations. With¬ 
out relying in dynamical simulations, Lyutikov, Pariev k 
Gabuzda (2005) considered the polarization properties of 
the synchrotron radiation emitted by relativistic, force-free 
jets with large scale helical fields and prescribed emissivity 
and internal structure, and Broderick k Loeb (2009) stud¬ 
ied the signatures of relativistic helical motion in the rota¬ 
tion measures of parsec-scale jets as a way to extract the 
intrinsic polarization angles, basing on a very simple heli¬ 
cal jet model. Following the same procedure as in this last 
work, Broderick k McKinney (2010) analyzed the rotation 
measures associated with jets formed from rapidly rotating, 
accreting black holes through the self-consistent general rel¬ 
ativistic MHD simulations of McKinney k Blandford (2009) 
extrapolated to parsec scales. 

The main goal of this paper is to build equilibrium mod¬ 
els of relativistic magnetised, axisymmetric jets with rota¬ 
tion propagating through an homogeneous, unmagnetised 
ambient medium at rest, which can serve as initial mod¬ 
els for ensuing magnetohydrodynamical simulations of com¬ 
pact extragalactic jets. Both the jet matter and the ambi¬ 
ent medium are described by an ideal gas equation of state 
with constant adiabatic index. As equilibrium solutions, our 
results complement those of Zakamska, Begelman k Bland¬ 
ford (2008) (for self-similar axisymmetric relativistic, non¬ 
rotating jet models with a purely toroidal magnetic field and 
an ultrarelativistic equation of state, in a pressure declining 
atmosphere), Gourgouliatos et al. ( 2012 ) (for non-rotating, 
relativistic magnetised jets without surface currents), Bodo 
et al. (2013) (who analyzed the stability of cold relativis¬ 
tic magnetised cylindrical flows) and Lyubarsky (1999) (on 
force-free jets). Beyond axisymmetric solutions, Mizuno et 
al. (2009, 2012) have investigated the influence of jet ro¬ 
tation and shear on the development of the current-driven 
kink instability of force-free helical magnetic flows via three- 
dimensional relativistic magnetohydrodynamic simulations. 
In the second part of the paper, we use these equilibrium 
solutions as initial models for pilot time-dependent RMHD 
simulations in one and two spatial dimensions. 

The organization of the paper is as follows. In Sect. 2 
the basic assumptions of our model are presented. The pres¬ 
sure equation leading to the transversal equilibrium for both 
rotating and non-rotating jet models is derived in Sect. 3. 
The remaining jet functions needed to solve the equilibrium 
equation are also defined in this section. In Sect. 4, the re¬ 
sulting equilibrium models are represented in the averaged 
magnetic pitch angle/averaged magnetisation plane and the 
forbidden regions in this plane identified. Section 5.1 fo¬ 
cusses on the capability of our RMHD code to maintain 


steady, one-dimensional solutions along dynamically signif¬ 
icant timescales. In Sect. 5.2, we present the first RMHD 
simulations of steady, two-dimensional relativistic, magne¬ 
tised, rotating jets with shocks. In Sect. 6 , the analytical 
solutions analyzed in the first part of the paper and the nu¬ 
merical simulations are put into an astrophysical context. 
A brief summary of the paper along with the most relevant 
conclusions are presented in Sect. 7. Finally, a short descrip¬ 
tion of the numerical RMHD code used in the simulations 
and its validation can be found in the appendices. 


2 MODEL ASSUMPTIONS 


In this paper, we seek solutions of steady, relativistic, mag¬ 
netised axisymmetric jets propagating through an unmag¬ 
netised ambient medium at rest. Units are used in which 
the light speed (c), the ambient density (pa) and the radius 
of the jet {Rj) are set to unity. Both the jet and the ambi¬ 
ent medium plasmas are assumed to behave as a perfect gas 
with constant adiabatic index 7 = 4/3. 

In order to make the problem tractable, we adopt 
several simplifications. As said in the previous paragraph, 
the jets are assumed to be axisymmetric. Using cylindri¬ 
cal coordinates (referred to an orthonormal cylindrical basis 
{er, e(/>, e^}) in which the jets propagate along the z axis, 
axisymmetry implies that there is no dependence on the az¬ 
imuthal cylindrical coordinate, cj). In the first part of the 
paper, we shall further assume that the jet models have slab 
symmetry along the z axis. This means that the radial mag¬ 
netic field, ^ is zero. Finally, this symmetry condition to¬ 
gether with the assumed stationarity of the flow, forces the 
radial velocity, v'^, to be zero too. Hence the jet solutions are 
characterised by six functions, namely the density and the 
pressure, p{r) and p(r), respectively, and the two remaining 
components of the velocity, v^{r), and of the magnetic 

field, B^(r). The ambient medium is characterised by 

a constant pressure, pa (besides pa = va — a = 0 ? 

Bl = Bt = Bl = 0). 

Under these conditions, the equation of transversal 
equilibrium establishing the radial balance between the to¬ 
tal pressure gradient, the centrifugal force and the magnetic 
tension, allows to find the equilibrium profile of one of the 
variables in terms of the others. We shall fix the radial pro- 
hles of p, ^ B^ and B^, and solve for the profile of the 


gas pressure, p. 

, We use top-hat prohles for p, 

and B^ 

p{r) = 1 

Pj, 

1, 

0 ^ r ^ I 
r > I, 

(1) 

II 

f 

[ 0 , 

0 ^ r ^ I 
r > I, 

(2) 

B^ir) = 

/ 

1 0. 

, O^r^I 
r > I, 

(3) 


(where pj , Vj and Bj are constants) and more complex pro¬ 
files for the two remaining functions. In particular, for the 
azimuthal magnetic held in the laboratory frame we choose 


B^{r) 


1 + ’ 

0 , 


0 ^ r ^ I 
r > I. 


(4) 
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This function represents a toroidal magnetic field that grows 
linearly for r <C reaches a maximum at r = 

m? then decreases as 1/r for r ^ RBd>,m ^tnd is set 
equal to zero for r > 1. It is a smooth fit of the piecewise 
profile used by Lind et ah (1989) (see also Komissarov 1999a; 
Leismann et ah 2005) and corresponds to a uniform current 
density for radius r <C Rb^ declining up to r = 1, and a 
return current at the jet surface. 

For the rotation profile we consider three situations: 

(i) Models without rotation: u'^(r) = 0. 

(ii) Models with rigid rotation: 


v^(r) 


0 ^ r ^ 1 
0 , r > 1. 


( 5 ) 


(iii) Differentially rotating models with a smooth tran¬ 
sition between an inner jet core with rigid rotation and a 
Keplerian sheath^: 


v^(r) 


0 , r > 1. 


( 6 ) 


Note that we associate the rotation of the jet with a 
non-zero azimuthal flow velocity, and not with the rotation 
of the magnetic field lines (u^ — j which can 

be non-zero even in the case of = 0. 

Our assumptions exclude by construction the important 
class of force-free equilibrium solutions considered by, e.g., 
Lyubarsky (1999); Mizuno et ah (2009, 2012), in which the 
gas pressure and the matter inertia are negligible. In the 
context of the present study, the force-free solutions should 
be understood as complementary. On the other hand, the 
analysis performed in this paper can be applied to any ra¬ 
dial profiles of density, and axial and toroidal flow veloc¬ 
ity and magnetic field, in particular to those derived from 
the self-consistently magnetically launched RMHD models 
in Komissarov et al. (2007, 2009) or in Lyubarsky (2010) 
papers. However these solutions are essentially multidimen¬ 
sional and do not belong to the class of solutions discussed 
in this work. 


3 TRANSVERSAL EQUILIBRIUM 

Under the conditions established in the previous Section, the 
RMHD equations in (orthonormal) cylindrical coordinates 
(see, e.g., the Appendix A in Leismann et al. 2005) reduce 
to a single ordinary differential equation for the transversal 
equilibrium of the jet 

^ ph*W\v^f - 

dr r ■ ^ 

In this equation, p* and h* stand for the total pressure and 
the specific enthalpy including the contribution of the mag¬ 
netic field 

p*=p+Y 

h* = 1 + £ + p/p + 5^/p, (9) 

^ Strictly speaking, the sheath will have a Keplerian rotation 
profile only for r > (or < 1). 

© 2015 RAS, MNRAS 000, 1-19 


where p is the fluid pressure, p its density and £ its spe¬ 
cific internal energy, {p = t, r, 0, z) are the components 
of the 4-vector representing the magnetic field in the fluid 
rest frame and stands for where summation over 

repeated indices is assumed^, (i = r, 0, z) are the compo¬ 
nents of the fluid 3-velocity in the laboratory frame, which 
are related to the flow Lorentz factor, W, according to: 


W = 




( 10 ) 


The following relations hold between the components 
of the magnetic field 4-vector in the comoving frame and 
the three vector components B^ measured in the laboratory 
frame: 


b° = WB%i , 


( 11 ) 

( 12 ) 


The square of the modulus of the magnetic field can be writ¬ 
ten as 


d 2 

j2 _ I / \2 

^ = vr + (^ 


(13) 


with B^ = B^Bi. 

Equation (7) establishes the transversal equilibrium be¬ 
tween the total pressure gradient and the centrifugal force 
(first term on the r.h.s.), that tends to produce a positive 
gradient of the radial total pressure profile, and the magnetic 
tension (second term on the r.h.s), that in turn favours to 
increase the total pressure towards the axis. Once fixed the 
radial profiles of p, , B"^ and B^, we solve the resulting 

first-order, linear, non-homogeneous equation for the profile 
of the gas pressure, p, together with the boundary condition 
at r = 1 given by pi — Pa^ where pi stands for p*(r= 1)^. 

Unless otherwise established, all the models considered 
in this paper have been computed for pj =0.01, = 0.97, 

Rb^ ~ 0.37, Pa — 0.1. 


3.1 Jets without rotation 


In the case of a jet without rotation, the equilibrium equa¬ 
tion (7) can be easily transformed into 


(^ 0)2 Qcf 


(14) 


dr rW^ VF2 dr ’ 

which can be integrated by separation of variables to give 

\ 2 


p{r) = < 


Bf 


W(1 + (r/i?B..^)2) 


1 . Pa, 


+ (7, 0 ^ r ^ 1 


r > 1. 


(15) 


Using the boundary condition, the integration constant C 
can be fixed to be 

?^^2 / d 0\2 


C ^ Pa- 


(Bt) 


^(i+(fls..j^). 


(16) 


^ Note that we are absorbing the (constant) magnetic permeabil¬ 
ity, 47r, in the definition of the magnetic field. 

^ Along this work, quantities with subindex 1 refer to their values 
at r = l, i.e., qi := q(r = l). 
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It can be easily seen that the presence of a toroidal 
field like the one defined in (4) increases the gas pressure 
up to r = Rb*+ 1) (~ 0-53, for 
Rb<i> m — 0-37) and decreases it outside so that the averaged 
gas pressure inside the jet remains unchanged with respect 
to the case of zero toroidal magnetic field (see next Section). 

Figure 1 displays two representative equilibrium models 
of non-rotating jets. In the left panel, the toroidal magnetic 
field is small enough as to produce an almost constant gas 
pressure profile inside the jet. In the case of the model dis¬ 
played in the right panel, the magnetic tension makes the 
gas pressure drop three orders of magnitude across the jet^ 
producing an equilibrium model with a central spine with 
high pressure (see inset panel). 


mean values of the toroidal magnetic field, the rotation ve¬ 
locity and the gas pressure, namely and according 

to 


_ fo q{r)rdr 


Qj ■- 


= 2 / q(r)rdr. 

Jo 


fo 

For the averaged toroidal magnetic field, we get 


^rctan 


(19) 


( 20 ) 


(ft! 0.81 for = 0.37). 

— 2 

In the case of models with rigid rotation, 

j 3 

whereas for the differentially rotating models. 


3.2 Jets with rotation 


For jet models with rotation, the equation for the transversal 
equilibrium (7) reads 


dp 

dr 


ct>\2 


r{W^ 


^ ^ ) dr 

^,dv^ , phW^v^f 


^ ^ dr 

(sv -2 bv) 


(17) 


where = (l — (r’^)^) and h — 1 + £ + p/p. Equa¬ 

tion (17) is solved numerically, with a standard fourth-order 
Runge-Kutta method, starting from the jet surface, where 


Pi =Pa 



{Bff 

2{w^y 


(18) 


Two illustrative equilibrium models of differentially ro¬ 
tating jets (^^m ~ — 0-25) are shown in Fig. 2. 

In the left panel, the gas pressure profile of the equilibrium 
model displays a deep minimum at r ~ 0.25 due to the cen¬ 
trifugal force caused by the jet rotation. Models with larger 
azimuthal speeds and/or smaller azimuthal magnetic fields 
would make the minimum pressure to reach negative val¬ 
ues. In the case of the model shown in the right panel, the 
toroidal magnetic field in the inner region is large enough to 
keep the pressure high, despite the action of the centrifugal 
force. 


4 EQUILIBRIUM MODELS IN THE 0-^ PLANE 

The equilibrium jet models can be displayed in a magnetic 
pitch angle/magnetisation (0 — /d) plane in terms of the 
toroidal and axial magnetic field components, for fixed jet 
density and kinematics (axial velocity and rotation velocity 
profile) and fixed ambient pressure. 

In order to compute the averaged magnetisation and 
magnetic pitch angle across the jet, we first compute the 


^ For the chosen dependence of with radius, Eq. (4), the gra¬ 
dient of the magnetic pressure is proportional to (and smaller 
than) the magnetic tension, with the proportionality factor rang¬ 
ing between (1 -h (1/R^0 (0.12 for ^ = 0.37) and 

1.0 across the jet radius. 



(ft 0.82for = 0.25). 

Einally, in the absence of rotation^, 

Pj =Pa - 


( 21 ) 


( 22 ) 


In the rotating cases, the averaged pressure in the jet is 
computed numerically. 

Now, the averaged pitch angle of the magnetic field, 
and the averaged magnetisation, /3j, are defined as 




( 


4^3 

:= arctan 


(23) 





Wj 

_ 3 

^Pj 


(24) 

In 

this last expression. 



b4> 

_ 



■ 

+ (B|)^ (l - 

(25) 


4.1 Jets without rotation 


Eigure 3 displays non-rotating equilibrium jet models in the 
magnetic pitch angle/magnetisation plane for fixed jet den¬ 
sity (pj = 0.01) and ambient pressure {pa = 0.1), and two 
axial jet velocities {vj = 0.5, left panel; Vj = 0.97, right 
panel). Drawn are lines of constant with increas¬ 

ing from left to right along each line (the radius at which 
the toroidal magnetic field reaches its maximum, 
has been fixed to 0.37 for all the models). There is a max¬ 
imum axial magnetic field, Rym = ^2pa, corresponding to 
a purely axial magnetic field compatible with a positive gas 
pressure within the jet. The first line starting from the top 
corresponds to (<^ = 0.975). Then, the ax¬ 

ial magnetic field labeling each line decreases linearly from 




top to bottom down to . Along each line, increases 


from left to right up to the largest value compatible with 
a positive gas pressure at the jet surface, for the given Bj. 


^ Note that pj for the non-rotating case is independent of the 
toroidal magnetic field, as advanced in the previous Section. 
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Figure 1. Pressure equilibrium profiles for two representative jet models without rotation. Model parameters: pj = 0.01, Vj = 0.97, 
~ ^ 10“^ (left panel), 2.81 (right panel), ^ = 0.37, = 0.436 (left panel), 4.36 x 10“^ (right panel), pa = 0.1. 


By construction, the region beyond the line corresponding 
to the maximum axial magnetic field can not be filled with 
equilibrium jet models within the class of solutions sought 
in this work (forbidden region I). 

Note that the maximum axial magnetic field, de¬ 

fined in the previous paragraph, is independent of the jet 
Lorentz factor. Hence, the left and right panels of Fig. 3 
display equilibrium models covering the same range of ax¬ 
ial magnetic fields. On the other hand, increasing the jet 
Lorentz factor seems to stretch the lines towards larger val¬ 
ues of the magnetic pitch angle. This is a consequence of the 
role played by the Lorentz factor in the magnetic pressure, 
reducing the relative weight of the axial magnetic field with 
respect to the toroidal one for given magnetisation. Finally, 
note that the region beyond the right ends of these Bj = 
constant lines towards the top-right corner of the plot is also 
forbidden (forbidden region II), i.e., there are no equilibrium 
models (within our class of solutions) available in this region 
of the plane (for given jet density and axial speed, and fixed 
ambient pressure). 


4.2 Jets with rotation 

Equilibrium models corresponding to rotating jets are dis¬ 
played in Figs. 4 and 5. Figure 4 corresponds to models with 
rigid rotation whereas Fig. 5 displays equilibrium models 
with an inner jet core with rigid rotation shrouded by a Ke- 
plerian sheath (junction radius, equal to 0.25). The 


remaining parameters are as in Fig. 3 (pj = 0.01, Vj = 0.97, 
= 0.1) exception made of the maximum azimuthal veloc¬ 
ity, which was set to 0.1 (left panels in Figs. 4 and 5) 
and 0.2 (right panels). Again as in Fig. 3, drawn are lines of 
constant with increasing from left to right along 

each line. 

The maximum axial magnetic field corresponding to a 
purely axial magnetic field compatible with a positive gas 
pressure within the jet is now 


Birr. = 


2pa 


1 - {vtr 


(26) 


Again, the first line starting from the top corresponds to 
Bj'm — (^Bj rn = 0.975). Then, the axial magnetic field 
labeling each line decreases linearly from top to bottom up 
B^' 

to . Along each line, Bf increases from left to right up 
to the largest value compatible with a positive gas pressure 
at the jet surface, for the given Bj. Let us note that, un¬ 
like the case of non-rotating jet models, there are now two 
families of equilibrium models depending on the sign of the 
product v^(r)B^(r). Figures 4 and 5 show models of the 
family with positive sign in which the helices corresponding 
to the magnetic field and the fluid stream lines turn in the 
same direction. Although similar concerning the transversal 
structure of the equilibrium models, both families can have 
qualitatively different emission imprints. The forbidden re¬ 
gions I (for axial magnetic fields beyond the maximum) and 
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Figure 2. Pressure equilibrium profiles for two representative jet models with the same rotation profile (rigidly rotating inner jet core 
and a Keplerian sheath). Model parameters: pj = 0.01, ^ = 0.20, R^ 4 , ^ = 0.25, = 0.97, ^ = 0.151 (left panel), 0.443 (right 

panel), R^cj) ^ = 0.37, = 0.396 (left panel), 4.41 x 10“^ (right panel), pa = 0.1. 


II (beyond the right ends of the Bj = constant lines towards 
the top-right corner of the plot) are also present in the case 
of rotating jets. 

Let us now consider the left panels of Figs. 4 and 5. 
The differences among the lines drawn in these panels and 
those displayed in the right panel of Fig. 3 are small. This 
is a sign that the azimuthal velocity is small and does not 
affect much the structure of the models. By contrast, the 
lines displayed in the right panels of Figs. 4 and specially 
in those of Fig.5, occupy a different region in the pitch an¬ 
gle/magnetisation diagram. This is mainly due to the fact 
that for rotating jet models there is another forbidden region 
beyond the left ends of the lines and towards the bottom-left 
corner. In this region (forbidden region III), models of the 
kind analyzed in this work would have sections with nega¬ 
tive gas pressure for given azimuthal velocity profile due to 
the centrifugal force. The effect is larger in the case of mod¬ 
els with the maximum azimuthal velocity inside the jet for 
which the centrifugal force is larger (for the same as 

those shown in the right panel of Fig. 5. For the same reason, 
the presence of a non-zero azimuthal velocity reduces effec¬ 
tively the maximum axial magnetic field compatible with 
positive gas pressures at the jet surface, hence modifying 
the forbidden region I. Models corresponding to line 1 in 
the right panels of Figs. 4 and 5 have negative gas pressures 
due to an excess of centrifugal force. 


5 STEADY RELATIVISTIC MAGNETISED 
JETS 

The analytical solutions discussed in the previous Section 
are interesting on their own but limited. On one hand, these 
solutions are restricted to infinite planar-symmetric jets in 
pressure equilibrium (which, in particular, hampers the de¬ 
velopment of radial components of the flow velocity and the 
magnetic field). On the other hand, although extragalactic 
jets are stable to a certain extent, they are far from being 
steady. Hence the need to combine these analytical steady 
solutions, that could serve as initial conditions, with dynam¬ 
ical (and emission) simulations. This strategy has demon¬ 
strated its success in the context of purely relativistic jets 
and can be naturally extended to RMHD jets (see the In¬ 
troduction) . 

In Sect. 5.1 the analytical solutions described in Sect. 4 
are used to probe the capability of our RMHD code to 
maintain steady solutions along dynamically significant 
timescales. This is a crucial necessary step if one wants to 
use dynamical simulations as a numerical laboratory and 
keep control on the parameters of the simulated models. As 
a first application, in Sect. 5.2, the RMHD code will be used 
to compute the steady solutions corresponding to overpres¬ 
sured jets. 

The numerical RMHD code used in these simulations is 
a conservative, finite-volume code based on high-resolution 
shock-capturing techniques. Its characteristics, as well as the 


© 2015 RAS, MNRAS 000 , 1-19 
















The structure of steady, relativistic, magnetised jets with rotation 7 




Magnetic pitch angle [deg] 


Figure 3. Jet magnetisation versus magnetic pitch angle in terms of ^ and Bj for equilibrium models without rotation. Drawn are 

lines of constant SJ, with B^^ increasing from left to right along each line up to the maximum value compatible with a positive gas 
pressure at the jet surface (r = 1), for the given axial magnetic field. Model parameters: pj = 0.01, Vj = 0.50 (left panel), 0.97 (right 
panel), R^cf, ^ = 0.37, pa = 0.1. Line i (= 1,..., 10) corresponds to models with B^ = 0.975S?^(11 — 'i)/10, where Bj^ = 0.447 is 
the maximum axial magnetic field compatible with a positive gas pressure within the jet (and zero toroidal magnetic field). Forbidden 
regions (see text for definitions) are also indicated. 


specific algorithms used are briefly described in Appendix A. 
Appendix B tests the code performance by means of several 
ID and 2D standard problems. 

5.1 Models with pressure equilibrium: Analytical 
vs. dynamical solutions 

In this case, the jet models have slab symmetry along the z 
axis making the problem to find steady solutions one dimen¬ 
sional. We have used the RMHD code in the ID radial cylin¬ 
drical coordinate to test the code ability to keep steady the 
initial equilibrium state. Two revealing cases have been cho¬ 
sen and their evolution shown in Figs. 6 and 7. The first case 
(Fig. 6 ) corresponds to a model with maximal axial magnetic 
field, differential rotation (rigidly rotating inner jet core and 
Keplerian sheath) and minimal toroidal magnetic field. The 
different panels in the figure display the profiles of gas (and 
total) pressure and density in logarithmic scale, and the 
three components of the fluid velocity and magnetic field. 
The initial equilibrium profile (thin filled line) and six pro¬ 
files corresponding to times t = 100,110,120,130,140,150® 
are displayed in each of the panels. The code keeps the initial 

® Due to relativistic aberration (the jet is fiowing along the axis 
at a speed of 0.97), the transverse signal crossing time of the jet 


steady solution (despite there is a small radial velocity of the 
order of 10 “® at both the position of the toroidal velocity 
maximum and at the jet/ambient medium transition) since 
the profiles at t 7 ^ 0 appear almost superimposed. Due to the 
slab symmetry of the problem (that the code keeps exactly), 
the radial component of the magnetic field is identically zero. 
The mean relative errors in the jet of the remaining quan¬ 
tities, ignoring the cell besides the axis (where the relative 
errors of the azimuthal components of the velocity and the 
magnetic field are < 3% and < 5%, respectively) and those 
at the jet/ambient medium transition (where the errors can 
be arbitrarily large), are displayed in Table 1. These errors, 
invisible at the plot scale, can be reduced by increasing the 
numerical resolution (200 numerical cells per jet radius in 
this case). The difference between the total energy per unit 
jet length in the initial analytical equilibrium model and the 
model at t = 150 is of 1.43%. 

The model shown in Fig. 7 is similar to that shown in 
Fig. 6 with the same (maximal) axial magnetic field but 
with rigid rotation and different (minimal) toroidal mag¬ 
netic field. The panel distribution and the lines plotted are 
the same as in Fig. 6 . The maximum of the azimuthal ve- 

radius is 4.1 time units. Hence the models shown in the figures 
correspond to tens of transversal light-crossing times. 
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Figure 4. Jet magnetisation versus magnetic pitch angle in terms of ^ and B? for equilibrium models with rigid rotation in which 
the helices corresponding to the magnetic field and the fiuid stream lines turn in the same direction. Drawn are lines of constant ^ 

with B'^^ increasing from left to right along each line up to the maximum value compatible with a positive gas pressure at the jet 

surface (r = 1), for the given axial magnetic field and rotation speed. Model parameters: pj = 0.01, v'j ^ = 0.10 (left panel), 0.20 
(right panel), Vj = 0.97, ^ = 0.37, pa = 0.1. Line z (= 1,..., 10) corresponds to models with Bj = 0.975SJ^(11 — 'i)/10, where 

Bj^ = 0.449 (left panel), 0.462 (right panel) is the maximum axial magnetic field compatible with a positive gas pressure within the 
jet for the given rotation speed (and zero toroidal magnetic field). Models corresponding to line 1 in the right panel have negative gas 
pressures due to an excess of centrifugal force and are not shown. Forbidden regions (see text for definitions) are also indicated. 


Table 1. Mean relative errors in the jet for the equilibrium models displayed in Figs. 6 (differentially rotating model, DR) and 7 
(rigidly rotating model, RR) at t = 150. 


Model 

P 

P 

v<t> 

yZ 

b4> 

B^ 

DR 

< 10-3 

< 2 X 10-^ 

< 4 X 10“® 

< 5 X 10-'^ 

< 10-® 

< 10-5 

RR 

< 3 X 10-^ 

< 2 X 10“"^ 

< 10-^ 

< 10-® 

< 5 X 10-5 

< 2 X 10-^ 


locity is now at the jet surface and the minimum in pressure 
is shifted towards larger radius. The radial velocity is again 
of the order of 10“^ at large times proving that the code 
has found a steady solution. The mean relative errors in 
the jet of the remaining quantities (the radial component 
of the magnetic field is zero), excluding the cell besides the 
axis (where the errors of the azimuthal components of the 
velocity and the magnetic field are < 7% and < 8%, respec¬ 
tively) and those at the jet/ambient medium transition, are 
displayed in Table 1. The difference between the total en¬ 
ergy per unit jet length in the initial analytical equilibrium 
model and the numerical one at t = 150 is now of 1.18%. 


5.2 Overpressured jet models: conical shocks and 
Mach discs 

In this Section, the equilibrium profiles discussed in Sect. 4 
are used as a boundary condition to inject the jets into a two 
dimensional domain representing an ambient medium with 
a pressure mismatch. In their attempt to reach again the 
equilibrium, the jets undergo sideways motions generating 
radial components of the flow velocity and the magnetic field 
that break the slab symmetry of the original jet model along 
the z axis. 

Four simulations of overpressured, rotating jets are pre¬ 
sented in this Section. The jets are injected through a 
nozzle of radius equal to 1 into an axisymmetric cylindri¬ 
cal domain with (r, z) G [0,4] x [0,30]. The evolution of 
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Figure 5. Jet magnetisation versus magnetic pitch angle in terms of and B? for equilibrium models with differential rotation in 

which the helices corresponding to the magnetic field and the fiuid stream lines turn in the same direction. Drawn are lines of constant 
with B^^ increasing from left to right along each line up to the maximum value compatible with a positive gas pressure at the jet 

surface (r = 1), for the given axial magnetic field and rotation speed. Model parameters: pj = 0.01, ^ = 0.10 (left panel), 0.20 (right 

panel), = 0.25, = 0.97, — 0-37, pa = 0.1. Line i (= 1,..., 10) corresponds to models with = 0.975B|^(11 — 2)/10, 

where B^ ^ = 0.448 (left panel), 0.452 (right panel) is the maximum axial magnetic field compatible with a positive gas pressure within 
the jet for the given rotation speed (and zero toroidal magnetic field). Models corresponding to line 1 in the right panel have negative 
gas pressures due to an excess of centrifugal force and are not shown. Forbidden regions (see text for definitions) are also indicated. 


Table 2. Parameters defining the overpressured models and corresponding averaged jet values. 


Model 

Pa 

K 

Pj 

Rotation 

law“ 


V. 

"1 




cj) 

^j 

j 

Pj 

<pj 

[deg] 

hi 

J012 

0.2 

2.5 

0.01 

RR 


0.15 

0.97 

0.37 

0.5 

0.5 

0.100 

0.407 

0.246 

39.14 

0.602 

J013 

0.2 

2.5 

0.01 

RR 


0.24 

0.97 

0.37 

0.5 

0.5 

0.160 

0.407 

0.033 

39.14 

4.83 

J085 

0.2 

5.0 

0.01 

DR 

0.25 

0.34 

0.92 

0.37 

0.5 

0.5 

0.280 

0.407 

0.358 

39.14 

0.503 

J086 

0.2 

5.0 

0.01 

DR 

0.25 

0.15 

0.92 

0.37 

0.5 

0.5 

0.123 

0.407 

0.741 

39.14 

0.214 


^ DR: differential rotation as defined in Eq. (6); RR: Rigid rotation. 


the flow in the domain is simulated with the RMHD code 
(see the Appendix) in 2D radial, axial cylindrical coordi¬ 
nates with a resolution of 80 (40) cells per jet radius in 
the radial (axial) direction. In order to disturb the ambi¬ 
ent medium as little as possible along the simulation, the 
domain (r, z) G [0,1] x [0,30] is initially filled with the an¬ 
alytical, injection solution. Reflecting boundary conditions 
are set along the axis (r = 0) and at the jet base outside the 


injection nozzle (r > 1, z = 0). Zero gradient conditions are 
set in the remaining boundaries. 

The basic parameters defining the models together with 
the corresponding averaged jet values are displayed in Ta¬ 
ble 2. The models, set up to be in equilibrium with an ambi¬ 
ent pressure are injected into an atmosphere with pres¬ 
sure Pa = Pa! where K is the jet overpressure factor. 
Models JO 12 and JO 13 are injected with a moderate over- 
presure factor of iF = 2.5 and rotate rigidly with maximum 
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Figure 7. Time evolution of an equilibrium model similar to that shown in Fig. 6 but with rigid rotation (and different toroidal magnetic 
field: ^ = 4.00 X 10“^). Drawn are the initial equilibrium profile (thin filled line) and six lines at t = 100 (thin dotted line), 110 

(thin dashed line), 120 (thin dotted-dashed line), 130 (thick filled line), 140 (thick dotted line) and 150 (thick dashed line). The radial 
computational domain spans the interval r G [0,4]. Numerical resolution: 200 cells per jet radius. 
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speed = 0.15 (Model J012) and 0.24 (J013). The rest 
of the parameters are the same in the two models. However, 
the different rotation profiles at injection generate different 
equilibrium profiles of the total pressure which are on the 
basis of the different flow structure of both models. Figure 8 
shows the distributions of several quantities at time t = 100 . 
In both cases the new equilibrium states are set through a 
series of conical fast-magnetosonic shocks however, whereas 
the separation between shocks in model JO 12 is about 15 
initial jet radii, in the case of model JO 13 the separation 
is much shorter (of only 4 initial jet radii). Related to this 
dissimilarity is also the change in the jet radius, which in the 
case of model JO 13 remains very close to 1 , but in the case 
of JO 12 oscillates between 1 and 1.5. This different structure 
can also point to differences in their stability properties in 
response to small perturbations. 

Models J085 and J086 have an overpressure factor of 
K — 5.0. In this case, the new equilibrium conditions can 
not be settled through a series of conical shocks and two 
Mach discs are formed instead. A snapshot of the two mod¬ 
els at t = 40 (soon before model J085 is completely decol- 
limated and disrupted) is shown in Fig. 9. Besides the dif¬ 
ference in the overpressure factor between the two sets of 
models (the relevant one to produce the switch between con¬ 
ical shocks and Mach discs), models J085 and J086 follow 
a differentially rotating law, with maximum rotation speed 
— 0-34 (Model J085) and 0.15 (J086). The effective¬ 
ness of Mach discs in decelerating the flow is clearly seen 
in the Lorentz factor panel of Fig. 9. However, whereas in 
the case of model J086, the flow accross the oblique magne- 
tosonic shock recollimates downstream and survives at least 
up to t = 100, model J085 is completely decollimated and 
disrupted. 


6 DISCUSSION 

The gist of the paper is the construction of models of planar 
(i.e., infinite), axisymmetric, rotating relativistic magnetised 
jets in transversal equilibrium made in Sect. 3. By solving 
the equation of transversal equilibrium, we have analyzed 
the roles of the magnetic pressure gradient, the magnetic 
tension and the centrifugal force in determining the profile 
of the gas pressure across the jet. Although our analysis has 
been performed for a particular choice of radial dependences 
of the jet density, velocity and magnetic field, it can be easily 
extended to more complex/realistic profiles. However, we 
should note that the region of the parameter space explored 
in our study leaves outside the class of force-free models 
(see, e.g., Lyubarsky 1999; Mizuno et al. 2009, 2012), which 
in the context of the present study, should be understood as 
complementary. 

An interesting result of our analysis is the existence of 
forbidden regions in the space of parameters where models 
of the class considered in our study could not be settled. 
These forbidden regions, represented in the magnetic pitch 
angle/magnetisation diagram in Sect. 4 are associated with 
the existence of maximum axial and toroidal magnetic field 
components compatible with the prescribed equilibrium con¬ 
dition at the jet surface, and/or an excess of centrifugal force 
producing gaps with negative pressures in the jet. Based 
on these general grounds, the existence of similar forbid¬ 


den regions can be expected for other radial dependencies 
within the class of planar, axisymmetric, magnetised jets 
in transversal equilibrium. Whether these forbidden regions 
can be filled at least partially with mixed MHD/force-free 
models (with proper radial profiles) needs further investiga¬ 
tion. 

One of the limitations of our approach is the assump¬ 
tion of the slab symmetry along the jet axis that excludes 
of our analysis all the solutions focussing in the collimation 
and acceleration of the jet and that have a dependence in the 
z direction (e.g., Komissarov et al. 2007, 2009; Lyubarsky 
2010). Hence our analysis on the transversal structure of jets 
seems better suited for models beyond the acceleration and 
collimation region. In this context, the solutions discussed 
in this paper can serve as initial injection conditions for ax¬ 
isymmetric jet simulations. (Magneto)hydrodynamical sim¬ 
ulations of jets have revealed as a very succesful approach to 
interpret the varied phenomenology of these objetcs (Gomez 
et al. 1997; Komissarov & Falle 1997; Agudo et al. 2001; 
Aloy et al. 2003; Mimica et al. 2009; Nakamura, Garofalo & 
Meier 2010; Forth et al. 2011; Fromm et al. 2012; Nakamura 
& Meier 2014; Mizuno et al. 2015). In particular, very long 
baseline interferometric observations of jets often suggest the 
presence of quasi-steady features (Jorstad et al. 2005; Lis¬ 
ter et al. 2013), interpreted as recollimation shocks. Besides 
this, a strong recollimation shock could be behind the nature 
of the millimeter-wavelength radio core in blazars (Marscher 
et al. 2008, 2010 ) in which the interaction of new super lu¬ 
minal components with this shock is invoked to explain the 
production of 7 -ray flares in these sources. The same kind 
of structure although parsecs away from the radio core has 
been suggested to explain the quasi-stationary bright radio 
feature in the jet of BL Lac located 0.26 mas from the core 
(Gohen et al. 2014, 2015), the component G80 in the 3G120 
jet (Roca-Sogorb et al. 2010) and the HST-1 complex in 
M87 (Giroletti et al. 2012 ). 

The simulations presented in Sect. 5.2 show steady jets 
with recollimation shocks of different characteristics. These 
simulations must be considered as exploratory since their 
purpose is only to probe the role of the overpressure fac¬ 
tor and the rotation speed of the jet in the properties of the 
internal magnetosonic shocks. According to this, it is impor¬ 
tant to note that the parameters used in these simulations 
were chosen to produce a variety of internal structures from 
periodic recollimation shocks of different strengths and spac¬ 
ing, to isolated Mach discs. The systematic study of struc¬ 
tures like those suggested in the previous paragraph and 
their observational signatures with a more physically ori¬ 
ented selection of parameters is one of the main objectives 
of our future work along this line of research. 


7 SUMMARY AND CONCLUSIONS 

In this paper, we present equilibrium models of relativis¬ 
tic magnetised, axisymmetric jets with rotation propagating 
through an homogeneous, unmagnetised ambient medium at 
rest. Under these conditions, the jet models are characterised 
by six functions defining the transverse profiles across the jet 
radius of the density and the pressure, p{r) and p(r), respec¬ 
tively, and the toroidal and axial components of the velocity, 
v^{r), r’^(r), and of the magnetic field, B^{r), B^{r) (the 
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Figure 8. Bottom halts of panels: Model J012. Top halts of panels: Model J013. Models parameters: see Table 2. The two contours 
in the top and bottom halts of the density map establish the transition between the jet (inner contour: 90% of jet material) and the 
ambient medium (outer contour: 10% of jet material). The steady states are set through a series of conical fast-magnetosonic shocks 
(with a spacing of about 15 initial jet radii in model J012, and of about 4 initial jet radii in model J013). 


corresponding radial components are zero by the impossed 
axisymmetry), and the ambient pressure, pa- In order to 
reduce the number of degrees of freedom of our study, we 
fix constant the jet rest-mass density and the axial compo¬ 
nents of the flow velocity and the magnetic field and ana¬ 
lyze the influence of a toroidal magnetic field and several 


rotation laws on the structure of the equilibrium models. 
The analysis performed in this paper can be applied to any 
radial profiles of density, and axial and toroidal flow veloc¬ 
ity and magnetic field, in particular to those derived from 
the self-consistently magnetically launched RMHD models 
of Komissarov et al. (2007, 2009) or Lyubarsky (2010), with 
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Figure 9. Bottom halts of panels: Model J085. Top halts of panels: Model J086. Models parameters: see Table 2. The two contours 
in the top and bottom halts of the density map establish the transition between the jet (inner contour: 90% of jet material) and the 
ambient medium (outer contour: 10% of jet material). Two Mach discs are clearly seen dX z ^ 7.0 (model J085) and 2 ; 12.5 (model 

J086). 
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qualitatively similar conclusions. However these solutions 
are essentially multidimensional and do not belong to the 
class of solutions discussed in this work. It is important to 
note that our assumptions exclude by construction the im¬ 
portant class of force-free equilibrium solutions considered 
by, e.g., Lyubarsky (1999); Mizuno et al. (2009, 2012), in 
which the gas pressure and the matter inertia are negligible. 
In the context of the present study, the force-free solutions 
should be understood as complementary. 

In the case of models without rotation (Sect. 3.1), the 
otherwise constant gas pressure profile is modified by the 
presence of the toroidal magnetic field that increases the 
pressure in a central spine around the jet axis due to the 
magnetic tension. The centrifugal force caused by a rotation 
velocity (Sect. 3.2) can change this structure producing a de¬ 
pression inside the jet. Models with large enough azimuthal 
speeds and/or small enough azimuthal magnetic fields would 
make the minimum pressure to reach negative values. 

In Sect. 4, the set of equilibrium models is represented in 
a magnetic pitch angle/magnetisation (0-yd) plane in terms 
of the toroidal and axial magnetic field components, for fixed 
jet density and kinematics (axial velocity and rotation ve¬ 
locity profile) and fixed ambient pressure. Several forbidden 
regions are identified in this diagram. First of all, there is 
a maximum axial magnetic field compatible with a positive 
gas pressure within the jet. Then, for given axial magnetic 
field, there is also a maximum toroidal magnetic field com¬ 
patible with a positive gas pressure at the jet surface. As 
a result, there is a forbidden region in the 0-yd diagram to¬ 
wards large magnetic pitch angles, large magnetisations. For 
rotating jet models there is an additional forbidden region 
towards the small pitch angle, small magnetisation corner 
in the 0-yd plane. Models in this region would have sections 
with negative gas pressure for given azimuthal velocity pro¬ 
file due to the centrifugal force. The effect is larger in the 
case of models with the maximum azimuthal velocity in¬ 
side the jet. The present study can be easily extended to jet 
models with different transversal profiles and magnetic field 
configurations. 

In Sect. 5 we have tested the ability of our RMHD 
code to maintain steady equilibrium models of axisymmet- 
ric RMHD jets in one and two spatial dimensions. The one 
dimensional simulations presented in Sect. 5.1 can be also 
understood as a consistency proof of the fidelity of the ana¬ 
lytical steady solutions discussed in this paper. Finally, the 
present study allows us to build equilibrium jet models with 
selected properties that could serve as initial conditions for 
dynamical (and emission) simulations of magnetised rela¬ 
tivistic jets with rotation, which will be the subject of fu¬ 
ture research. In particular, the exploratory simulations pre¬ 
sented in Sect. 5.2 were designed to probe the role of the 
overpressure factor and the rotation speed of the jet in the 
properties of the internal recollimation magnetosonic shocks, 
invoked to explain several observational trends in parsec- 
scale extragalactic jets. 
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i) Gell reconstruction: second-order accurate values of the 
primitive variables V = (p,p, at the left and right 

ends of the cells are obtained with linear functions and sev¬ 


eral limiters (MINMOD, VAN LEER, MG). MG and VAN 
LEER limiters can be degraded to MINMOD in case of 
strong shocks. No jump is allowed in the normal component 
of B at a cell boundary and the corresponding staggered 
magnetic field is used. 

ii) Riemann solvers: intercell numerical fluxes are com¬ 
puted by means of HLL and HLLG (Mignone V Bodo 2006) 
Riemann solvers. Accurate bounds of the maximum speeds 
of left and right propagating waves are obtained by solving 
the corresponding characteristic equation for the left and 
right states of each numerical interface. 

iii) Time advance: the multidimensional equations of 
RMHD are advanced in time in an unsplit manner using 
TVD-preserving Runge-Kutta methods of second and third 
order (Shu V Osher 1988, 1989). The time step is determined 


according to At = CFL x min 



(ID, cylindrical ra- 


1-1 . N A. CFL . f Ar Az \ 

dial version) or to At = —— x mm —-r, —-r (2D, 

V2 V|Ar,z,i| \^z,i,j\ J 

cylindrical axisymmetric version), where Xr,ij and Xz,i,j are 
the speeds of the fastest waves propagating in cell z, j along 
the r and z direction, respectively. 

iv) Gonstrained transport scheme as in Balsara V Spicer 
(1999) for multidimensional calculations. The conserved 
variables can be corrected with the algorithms described in 
Marti (2015). 

v) Primitive variables are recovered as in the IDvc method 
of Noble et al. (2006) and solving the resulting equation in 
Z = phW^ by bisection. 


The code advances the total energy density without the 
rest-mass energy density. This strategy improves the perfor¬ 
mance of the conservative scheme when the total energy is 
dominated by the rest-mass energy. 

The ID calculations presented in Sect. 5.1 have been 
performed with the MG limiter, and the HLLG Riemann 
solver, whereas those shown in Sect. 5.2 have been done with 
the VAN LEER limiter degraded to MINMOD at shocks, 
and the HLL Riemann solver. The advance in time has been 
done using the third-order TVD-preserving Runge-Kutta 
with CFL — 0.6 in all the cases. The relativistic correction 
algorithm GA2’ of Marti (2015) has been used to correct the 
conserved variables after each time step in the 2D simula¬ 
tions. 


APPENDIX B: CODE VALIDATION 

The code passes all the standard RMHD tests in ID (Komis¬ 
sarov 1999a; Balsara 2001; Giacomazzo V Rezzolla 2006) 
and 2D with the expected accuracy and convergence rate. A 
complete testing of the code will be presented elsewhere. 


APPENDIX A: RMHD CODE SUMMARY 

Our RMHD code is a conservative, second-order, finite- 
volume, constrained-transport code based on high- 
resolution, shock-capturing techniques. Its basic ingredients 
are the following: 


B1 Code accuracy 

The nominal second order of accuracy of the code is verified 
by means of the ID smooth test proposed by Del Zanna et 
al. (2007) describing the propagation of a large-amplitude, 
circularly polarized Alfven wave along a uniform background 
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field Bq. Taking the background field along the x axis, the 
transverse velocity components are 

Ott 27r 

v'^ = — AcOs[-—(x — Vat)], = — Asin[-—(x — Vat)] (Bl) 
A A 

(where A is the amplitude of the wave, and A its wavelength), 
and 


= -Bov^va, = -Bov^va. 


(B2) 


In the previous expressions, the speed of the Alfven wave, 
Va, is given by: 


Va = ^ 


gg(i - 

poho + B^{1 - A^)’ 


(B3) 


where po and ho are the density and the specific enthalpy of 
the background uniform medium^. 

Del Zanna et ah (2007) utilized this test to assess the or¬ 
der of accuracy of their code ECHO, and Beckwith & Stone 
(2011) did the same for the RMHD module of ATHENA. We 
have tested the accuracy of our code by measuring the errors 
on one of the transversal velocity components, namely v^, 
at one period, t = X/va, compared to the initial condition 
at t = 0. Our simulations cover the spatial domain [0,27r], 
corresponding to a wavelength, with periodic boundary con¬ 
ditions at both ends. Eor the background state we choose 
Po = 1, ho = 5 (7 = 4/3), Ho = 1- The amplitude of the 
wave is taken as 


A = 


7 + 3V5’ 


All the values, including the peculiar value of A, have been 
chosen to define the same wave as the one originally used by 
Del Zanna et al. (2007). 

Table B1 shows the errors and convergence orders in the 
Li norm (the absolute error averaged over the whole com¬ 
putational domain) for the test at various resolutions with 
MINMOD, VAN LEER and MC limiters. The HLLC Rie- 
mann solver and the third-order Runge-Kutta with CFL = 
0.6 has been used in all the cases. 


Table Bl. Accuracy of the code from the circularly polarized 
Alfven wave test. 


Method 

N 

Li 

error 

Li order 

MINMOD 

8 

1.83 

X 

10-1 

— 


16 

8.16 

X 

10-2 

1.17 


32 

2.19 

X 

10-2 

1.90 


64 

5.47 

X 

10-3 

2.00 


128 

1.51 

X 

10-3 

1.86 


256 

4.05 

X 

lO-i 

1.90 


512 

1.05 

X 

10-^ 

1.95 

VAN LEER 

8 

1.49 

X 

10-1 

— 


16 

3.96 

X 

10-2 

1.91 


32 

7.86 

X 

10-3 

2.33 


64 

1.58 

X 

10-3 

2.31 


128 

3.55 

X 

10-^ 

2.15 


256 

8.36 

X 

10-5 

2.09 


512 

2.01 

X 

10-5 

2.07 

MC 

8 

1.26 

X 

10-1 

— 


16 

2.70 

X 

10-2 

2.22 


32 

5.45 

X 

10-3 

2.31 


64 

1.28 

X 

10-3 

2.09 


128 

3.13 

X 

10-4 

2.03 


256 

7.75 

X 

10-5 

2.01 


512 

1.93 

X 

10-5 

2.01 


and the contact discontinuity. Eigure Bl shows the numeri¬ 
cal results on this test at t = 0.4 using 1600 numerical cells 
in the domain x G [—0.5, 0.5], with the MC limiter, HLLC 
Riemann solver and Runge-Kutta of third order for time ad¬ 
vance. The numerical solution is stable and the high-density 
shell propagating to the right is captured quite accurately 
with the present resolution. The largest errors are found in 
the shell values of the transverse components of the flow ve¬ 
locity and, specially, of the magnetic field. Overall our results 
are only improved by those obtained by Mignone V Bodo 
(2006), who used a second-order, MUSCL-Hancock scheme. 


B2 Simulation of discontinuous solutions and thin 
structures 

The capability of our code in simulating flows involving dis¬ 
continuities is demonstrated by means of the Riemann prob¬ 
lem number 3 proposed by Balsara (2001), with initial data 
(p,p, x", x", H", B^)l = (1, 1000 , 0 , 0 , 0 , 10 , 7, 7), 

(p,p,x"^,x^,x^H"^,H^,H^)i^ = (1,0.1,0,0,0,10,0.7,0.7) at 
both sides of a discontinuity at x = 0. The adiabatic in¬ 
dex of the ideal gas equation of state is 7 = 5/3. Accord¬ 
ing to the analytical solution computed by Giacomazzo V 
Rezzolla (2006), the test develops a pair (slow and fast) 
of left-propagating rarefaction waves, and a pair of right- 
propagating shocks, separated by a contact discontinuity®. 
A high-density shell forms between the leading slow and fast 
right-propagating shocks (moving at almost the same speed) 

^ Note that that expression for Va is different from the one in 
Eq. (85) of Del Zanna et al. (2007), though equivalent. 

® As in the rest of coplanar Riemann problems, no Alfven waves 
develop after the break-up of the initial discontinuity. 


B3 Multidimensional problems 

Results on the cylindrical magnetised blast wave test 
(Komissarov 1999a), with magnetisations larger than 10®, 
can be found in Marti (2015). In order to test the code in 
cylindrical coordinates, we have considered a similar prob¬ 
lem simulating the blast of a spherical region embedded in 
a uniform, magnetised medium (see Mignone et al. (2007)). 
In this test, a sphere of gas with density p = 10~^ and pres¬ 
sure p = 1 , is embedded in a static uniform medium with 
p = 10“^ and p = 3 X 10“^. The sphere, with radius r = 0.8, 
is centered at the origin and a linear smoothing is applied 
for 0.8 ^ r ^ 1. The whole region is threaded by a constant 
vertical field in the ^-direction, = 1. The adiabatic index 
of the ideal gas equation of state is 7 = 4/3. 

We have run this test in cylindrical coordinates in the 
domain (r,z) G [0,6] x [— 6 , 6 ] with 512 (radial) x 1024 
(axial) numerical cells. Open boundary conditions are used 
along the boundaries of the computational doain exception 
made of the axis (r = 0 ), were reflecting boundary condi¬ 
tions are used instead. The same numerical ingredients of 
the code as for the 2 D jet simulations presented in Sect. 5.2 
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were chosen. Figure B2 shows the distributions of several 
representative quantities at t = 4.0. Our results can be 
directly compared with those presented by Mignone et al. 
(2007) in their Fig. 11. The difference in pressure between 
the spherical region and the ambient medium produces the 
expansion of the central region delimited by a fast magne- 
tosonic shock propagating radially at almost the speed of 
light. Because of the strong sideways magnetic confinement 
an elongated structure develops in the z direction with a 
maximum Lorentz factor of W ~ 4.6. This problem is par¬ 
ticularly challenging because of the very large magnetization 
/3 = 1.67 X 10^. The problem also serves to test the treat¬ 
ment of the geometrical source terms in the code and also 
the CT scheme in cylindrical coordinates. 
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Figure B2. Proper rest-mass density, gas pressure and magnetic pressure (in logarithmic scale), flow Lorentz factor and magnetic 
field components at t = 4.0 for the spherical magnetized blast wave test discussed in the text. The magnetic held lines are printed on 
top of the Lorentz factor plot, logp G [-5.82,-2.58], logp G [-4.59,-0.68], logprn ^ [-0.57,-0.16], W G [1.00,4.57], \B'^\ G [0,0.52], 
\B^\ G [0.74,1.19]. 
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